Objectifs de la formation

  1. Comprendre le séquençage Illumina et les fichiers FASTQ
  2. Manipuler des données de séquençage avec R
  3. Découvrir BARQUE, un pipeline de métabarcoding québécois
  4. Appliquer les concepts sur des exemples pratiques

Avant de commencer

Installation des librairies R requises

  1. Ouvrir RStudio
  2. Envoyer les commandes suivantes dans la console
# Installation des packages Bioconductor
if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install(c("ShortRead", "Biostrings", "Rqc"))

Introduction

Comment le séquencage est effectué?

Crédit: J. Couillard

DOI: 10.1016/B978-0-12-802234-4.00002-1

Comment le séquencage est effectué?

Note

  • Chaque nucléotide à un fluorophore spécifique (4 canaux).
  • Certains appareils (ex. NextSeq) vont utiliser deux canaux car c’est plus rapide et moins coûteux.

DOI: 10.1016/B978-0-12-802234-4.00002-1

Le signal renvoyé par le laser est traduit en chromatographe pour chaque amplicon séquencé.

Comment le séquencage est effectué?

Exemple de chromatographe

Source: LabXchange.org

Comment le séquencage est effectué?

Exemple de chromatographe

Important

Le score est basé sur:

  • La séparation entre les clusters (pureté du signal)
  • Le rapport signal/bruit
  • La distribution de l’intensités des 4 canaux
  • La position dans le read (la qualité diminue généralement vers la fin)

Source: wikipedia.org

Fichiers de sortie du séquenceur Illumina

Exemple de nomenclature de fichier de sortie du séquenceur Illumina: sample1_1_L001_R1_001.fastq.gz

Signification Exemple
sample1 Identifiant de l’échantillon sample1, sample2
_1 Numéro de réplicat 1, 2, 3…
L001 Numéro de ligne / lane: L001-L008
R1/R2 Direction de lecture R1=avant, R2=arrière
001 Segment de fichier 001, 002…

Note

  • Séquençage paired-end = fichiers R1 + R2 pour chaque échantillon avec l’ensemble des amplicons.
  • Numéro de Ligne ou lane = Une flowcell Illumina est divisée en plusieurs voies physiques parallèles où le séquençage se déroule simultanément.

Structure du fichier FASTQ

Contient plusieurs lignes dont un ensemble de 4 lignes correspond à un amplicon.

system.file(package="ShortRead", "extdata/E-MTAB-1147/ERR127302_1_subset.fastq.gz") |> 
  gzfile() |> 
  readLines(n = 12)
 [1] "@ERR127302.8493430 HWI-EAS350_0441:1:34:16191:2123#0/1"                  
 [2] "GTCTGCTGTATCTGTGTCGGCTGTCTCGCGGGACATGAAGTCAATGAAGGCCTGGAATGTCACTACCCCCAG"
 [3] "+"                                                                       
 [4] "HHHHHHHHHHHHHHHHHHHHEBDBB?B:BBGG<DDAA?AABFEFBDBD@DDECEE3>:?;@@@>?=BAB?##"
 [5] "@ERR127302.21406531 HWI-EAS350_0441:1:88:9330:2587#0/1"                  
 [6] "CTAGGGCAATCTTTGCAGCAATGAATGCCAATGGGTAGCCAGTGGCTTTTGAGGCCAGAGCAGACCTTCGGG"
 [7] "+"                                                                       
 [8] "IIIIHIIIGIIIIIIIHIIIIEGBGHIIIIHGIIHIIIIIIIHIIIHIIIIIGIIIEGIIGBGE@DDGGGIG"
 [9] "@ERR127302.22173106 HWI-EAS350_0441:1:91:10434:14757#0/1"                
[10] "TGGGCTGTTCCTTCTCACTGTGGCCTGACTAAAACCCAGTGGCATTAAGAAAGAGTCACGTTTCCCAAGTCT"
[11] "+"                                                                       
[12] "GGHBHGBGGGHHHHDHHHHHHHHHFGHHHHHHHHHHHHHHHHHHHHHGHFHHHHHHHHHHHHHH8AGDGGG>"

Structure du fichier FASTQ

system.file(package="ShortRead", "extdata/E-MTAB-1147/ERR127302_1_subset.fastq.gz") |> 
  gzfile() |> 
  readLines(n = 4)
[1] "@ERR127302.8493430 HWI-EAS350_0441:1:34:16191:2123#0/1"                  
[2] "GTCTGCTGTATCTGTGTCGGCTGTCTCGCGGGACATGAAGTCAATGAAGGCCTGGAATGTCACTACCCCCAG"
[3] "+"                                                                       
[4] "HHHHHHHHHHHHHHHHHHHHEBDBB?B:BBGG<DDAA?AABFEFBDBD@DDECEE3>:?;@@@>?=BAB?##"
  1. Ligne [1] : Identifiant de séquence (commence par @)
  2. Ligne [2] : Séquence d’ADN (A, T, G, C) de l’amplicon
  3. Ligne [3] : Séparateur (+)
  4. Ligne [4] : Scores de qualité (un par base)

Ligne d’identification

Décomposition de l’identifiant

@SH00321:6:BWR98207-2813:1:1101:1065:1015 1:N:0:AACCATAGAA+GGCGAGATGG
  • SH00321 : ID de l’instrument
  • 6 : Numéro d’exécution
  • BWR98207-2813 : ID de la flowcell
  • 1:1101:1065:1015 : Tuile et coordonnées X:Y
  • 1:N:0 : Numéro de lecture, flag de filtre, bits de contrôle
  • AACCATAGAA+GGCGAGATGG : Codes-barres d’échantillon (double indexation)

Interprétation des scores de qualité

Exemple de chaine de caractères donnant le score de qualité (1 caractère par nucléotide): GGGGGGGGGGGGGGGGGGGGGGGG9GG-G9G9G9GGGG

Caractère ASCII Score Phred Taux d’erreur Précision
G 38 0,016% 99,984%
9 24 0,40% 99,60%
- 12 6,31% 93,69%

Tip

Règle générale : Q30 ou plus = haute qualité (99,9% de précision)

Scores de qualité du séquencage

Les scores de qualité utilisent l’encodage ASCII.

 Quality encoding: !"#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHI
                   |         |         |         |         |
    Quality score: 0........10........20........30........40   

Scores de qualité du séquencage

Avec R, comment convertir l’encodage ASCII en score de Phred et en probabilité d’obtenir une erreur (précision)?

# Convertir la chaîne de qualité en scores Phred
quality_string <- "GGGGGGG99GG-GGG"
quality_scores <- as.integer(charToRaw(quality_string)) - 33

# Afficher la correspondance
data.frame(
  Caractere = strsplit(quality_string, "")[[1]],
  Score_Phred = quality_scores,
  Precision = paste0(round((10^(-quality_scores/10)) * 100, 2), "%")
)

Scores de qualité du séquencage

   Caractere Score_Phred Precision
1          G          38     0.02%
2          G          38     0.02%
3          G          38     0.02%
4          G          38     0.02%
5          G          38     0.02%
6          G          38     0.02%
7          G          38     0.02%
8          9          24      0.4%
9          9          24      0.4%
10         G          38     0.02%
11         G          38     0.02%
12         -          12     6.31%
13         G          38     0.02%
14         G          38     0.02%
15         G          38     0.02%

Lire un fichier FASTQ avec R

# Lire le fichier FASTQ
fq <- system.file(package="ShortRead", "extdata/E-MTAB-1147/ERR127302_1_subset.fastq.gz") |> ShortRead::readFastq()

# Examiner la structure
fq
class: ShortReadQ
length: 20000 reads; width: 72 cycles

Extraire les informations principales

# Nombre d'amplicons
length(fq)
[1] 20000
# Extraire les amplicons
sequences <- ShortRead::sread(fq)

# Extraire les identifiants
ids <- ShortRead::id(fq)

# Extraire les scores de qualité
qualities <- Biostrings::quality(fq)

# Premier identifiant et séquence
ids[1]
BStringSet object of length 1:
    width seq
[1]    53 ERR127302.8493430 HWI-EAS350_0441:1:34:16191:2123#0/1
sequences[1]
DNAStringSet object of length 1:
    width seq
[1]    72 GTCTGCTGTATCTGTGTCGGCTGTCTCGCGGGAC...GTCAATGAAGGCCTGGAATGTCACTACCCCCAG

Résumé statistique d’un fichier FASTQ

# Longueur des séquences
seq_lengths <- Biostrings::width(sequences)
summary(seq_lengths)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
     72      72      72      72      72      72 
# Distribution des scores de qualité
qual_matrix <- as(qualities, "matrix")
mean_quality <- rowMeans(qual_matrix)
summary(mean_quality)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
   2.00   33.53   37.35   34.84   38.85   39.99 
# Composition en bases
Biostrings::alphabetFrequency(sequences, as.prob = TRUE)[1:4, 1:4]
             A         C         G         T
[1,] 0.1944444 0.2638889 0.2916667 0.2500000
[2,] 0.2361111 0.2222222 0.3194444 0.2222222
[3,] 0.2361111 0.2638889 0.2222222 0.2777778
[4,] 0.2500000 0.2638889 0.1666667 0.3194444

Assurance qualité avec la librairie Rqc

fq_files <- system.file(package="ShortRead", "extdata/E-MTAB-1147/") |> list.files(full.names = TRUE)
qa <- Rqc::rqcQA(fq_files)
Rqc::rqcCycleQualityPlot(qa)

Analyser la qualité du séquencage

Pourquoi la qualité des lectures décroit au fur et à mesure des itérations?

  1. Dégradation des fluorophores
  2. Accumulation d’erreurs de phase, soit la désynchronisation des clusters
  3. Épuisement des réactifs, compétition pour les nucléotides par exemple
  4. Dommages à l’ADN par stress chimique par exemple
  5. Accumulation de sous-produits (ex. Résidus de réaction s’accumulent)

Points clés à retenir

  • FASTQ = format standard pour les données de séquençage
  • 4 lignes par lecture d’amplicon : ID, lecture de l’amplicon, séparateur, qualité de la lecture
  • Paired-end : les fichiers R1 et R2 contiennent les lecture d’amplicons à apparier
  • Scores de qualité : Plus élevé = meilleur (viser Q30+)
  • La qualité se dégrade vers la fin des lectures (normal)

Qu’est ce que BARQUE?

Un large éventail d’outils en bioinformatique

Hakimzadeh et al, 2024

Qu’est ce que BARQUE?

  • Programme écrit par Éric Normandeau (IBIS, ULaval)
  • Disponible sur Github depuis 2017
  • C’est un pipeline de traitement de données de séquences
  • Il agit comme chef d’orchestre d’un ensemble de programme

Qu’est ce que BARQUE?

Creative Commons
  • FLASH v1.2.11+: fusion des lectures paired-end
  • VSEARCH v2.14.2+: outil d’analyse de séquences (version obligatoire)
  • Trimmomatic: outil de filtrage et découpage des lectures
  • NCBI BLAST: annotation des séquences non-correspondantes
  • Un ensemble de programmes UNIX comme cut, sort, uniq etc. qui manipulent des chaînes de caractères ou des nombres
  • R & Pyhton: Visualisation et synthèse des résultats

Les grandes étapes de BARQUE

  1. Vérification des fichiers des entrées dans le
  2. TRIM - Nettoyage/découpage des séquences
  3. MERGE - Fusion des lectures paired-end (R1 et R2)
  4. Split Amplicon - Séparation des amplicons (TODO ASV vs OTUs)
  5. Retrait des chimères - Détection et suppression des séquences chimériques
  6. Assignement taxonomique - Attribution des espèces

Télécharger BARQUE

Pratique

  1. On se

Les étapes en tant qu’utilisateur

  1. On place les fichiers fastq.gz avec la bonne nomenclature dans le dossier

Types de primer

Test

Espèces intégré à la base de données par PRIMER

COI

Aller chercher la base de données NCBI pour les primers COI

Pourquoi utiliser docker

BARQUE requière plusieurs programmes qui ne sont pas tous compatibles avec l’environnement Windows

Hakimzadeh et al, 2024

Comment docker fonctionne avec -v

Mettre un schéma

Étape 1. Télécharger l’image docker

Étape 2. Execute l’image